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Q '_ We introduce a general theoretical scheme for a class of phenomena char- 

acterized by an extremal dynamics and quenched disorder. The approach is 

>' 
^O ' based on a transformation of the quenched dynamics into a stochastic one 

-. ' with cognitive memory and on other concepts which permit a mathematical 

O 

t^^ ' characterization of the self-organized nature of the avalanche type dynamics. 

In addition it is possible to compute the relevant critical exponents directly 

Vh ' from the microscopic model. A specific application to Invasion Percolation is 

pH ■ presented but the approach can be easily extended to various other problems. 

o' 
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Extremal models of self-organized criticality (SOC 0) have recently attracted significant 
theoretical attention. This class of models describes phenomena ranging from fluid displace- 
ment in porous disordered media [jl| , to punctuated biological evolution . In these models, 
at each time step the dynamical activity is concentrated on the site with the extremal value 
of a quenched disordered variable. This rule leads to a rich and complex behaviour, which 
has been widely studied @. 

In this letter we describe a general theoretical approach which addresses the ba- 
sic problems of extremal models: (i) the understanding of the scale-invariance and self- 
organization; (ii) the origin of the avalanche dynamics and (iii) the computation of the rele- 
vant critical exponents. We apply it specifically to Invasion Percolation (IP), but it can be 
easily extended to other models of this type like the Bak and Sneppen model [|1^ . 



The usual real space methods, like Fixed Scale Transformation [0, cannot address di- 
rectly the problem of the irreversible dynamics with quenched disorder. In order to overcome 
this basic problem we introduced a mapping of a quenched extremal dynamics into a stochas- 
tic one with cognitive memory, the quenched- stochastic transformation, also called Run Time 
Statistics (RTS) 0]. This approach was improved in various steps p,|,0 and now we can 
develop it into a general theoretical scheme. Its essential points are: 

- Quenched-stochastic transformation. 

- Identification of the microscopic fixed point dynamics (SOC). 

- Identification of the scale invariant dynamics for block variables. 

- Definition of /oca/ growth rules for the extremal model. This clarifies the origin of avalanche 
dynamics. 

- Use of the above elements in the FST scheme to compute analytically the relevant expo- 
nents of the model. 

Let us start with the quenched-stochastic transformation. We will discuss it for the case 
of bond Invasion Percolation (IP) ||T|. The IP model describes the capillary displacement of 
a fluid in a porous medium. The medium is represented as a network of bonds. To each 
bond i is assigned a quenched random number Xi extracted from a given distribution with 
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density p{x) = 1. The invading cluster evolves by occupying the bond with the smallest Xj 
on its perimeter. The basic idea is to map the deterministic extremal IP dynamics into an 
annealed stochastic process. A general stochastic process is based on the following elements: 
a) a set of time dependent dynamical variables {^j,*}; b) a Growth Probability Distribution 
(GPD) for the single growth step {fii,t}, obtained from the {?7i,t}; c) a rule for the evolution 
of the dynamical variables rji^t —^ Vi,t+i- 

Therefore, in order to map IP onto a stochastic process we have to identify the above 
elements. 

An insight into the essence of the question is given by by the following example. Consider 
two independent random variables Xi, X2 uniformly distributed in [0, 1] and let us eliminate 
the smallest, for example X2. The probability that X2 < Xi is 1/2. At the second "time 
step", we compare the surviving variable Xi with a third, uniform, random variable X3 
just added to the game and, again, we eliminate the smallest one. In this case we need to 
compute the probability yUs that X3 < Xi given that X2 < Xi. This, using the rules of 
conditional probability, reads: 



/is = P(X3 < Xi) = P(X3 < X1IX2 < Xi) 

P(X3<XinX2<Xi) 2 



(1) 



P(X2<Xi) 3' 

where P{A\B) is the probability of the event A, given that B occurred, and P{A{^B) is 
the probability of occurrence of both A and B. Equation |l| tells us that the information 
X2 < Xi changes in a conditional way the effective probability density pi{x) of Xi. In fact, 
by imposing the condition X2 < x (given x < Xi < x + dx) we get pi{x) = 2x. Qualitatively, 
the event X2 < Xi decreases the probability that Xi has small values. 

The above example contains the essential idea of the quenched- stochastic transformation. 
Let us come back to the IP model. In view of the above example, each perimeter variable 
will have, at time t, an effective density depending on the past growth history of the bond. 
If a variable has lost many times, it will have a density more and more concentrated on 
great values. The past history of a variable can be represented by only one parameter, the 
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number of time steps that the variable has spent in the perimeter, that we call "age" k 
of the variable. Variables with the same age k will have the same effective density Pk,t{x) 
{Po,t{x) = p{x)) at time t. At this point we can express the probability that a variable of age 
k is selected at time t, i.e. it is the smallest perimeter variable, given the effective densities 
Pk,t{x) of all the perimeter variables for a realization of the process at time t, in the following 
way: 

l^k,t = f'dxpk,t{x)Y[{l - PeA^))""'''"'', (2) 

Jo g 

where Pk,t{x) = Jq dypk^y), the product is intended over all the ages of the variables and 
ngt is the number of active variables of age 9 at time t. The product inside the integral 
takes into account the competition of the selected variable with each one of the other active 
variables, while the integral between and 1 takes into account all the possible values of 
the growing variable. The temporal evolution of the densities of the still active variables is 
then given by: 

P9+i,t+i[x)=pe,t[x) , A y. (3) 

Jo 1 - Pe,t[y) 

where mk^y) = [pk,t{^) Ue (1 - Pe^tix))'""'''^'-"] /iJ,k,t- 

Equations 0, |^ describe a quenched extremal process as a stochastic process with memory. 
The presence of memory is enlighted by the dependence of the GPD on k. It has been 
shown by M. Marsili 0] that a mean field like expansion of Eq.(^ in the limit t -^ oo gives: 
f^k,oo ~ 7Ffi>" '^^^ power law behaviour oi fik,t guarantees that screening is preserved at all 
scales, which is the condition to generate fractal structures [p!0|] . 

An important quantity to study in IP is the histogram $i(x), which is the distribution 
at time t of the perimeter variables. It has been shown numerically |jl| that the histogram 
of IP self-organizes asymptotically into a stable shape, a theta function with a discontinuity 
aX X = Pc < 1, where Pc is the critical bond percolation threshold. In ref. [Q an equation for 
the temporal evolution of $i(a;) is derived analytically directly from the RTS microscopic 
dynamics, without a priori assumptions. The final form of this histogram equation, obtained 
after some technical approximations [Q is: 
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d^Mx) = pnt<^Kx) 



"* -M-] 



(4) 



Ut + 1 

where Nt is the number of perimeter variables a time t for a reahzation of the growth 
process, Qt = i^t), ^t = (^t+i — Nf), the mean (...) is over different reahzations, and /5 is 
the solution of :/5 = 1 — e~^^^ . The solution of eq. § is shown to converge asymptotically to 
a theta function (fig.|l|): 

lim $i(x) = e{x - Pc) (5) 

i^oo l-pc 

where pc = 1/2 for 2 — d bond IP. This result clarifies the SOC nature of the problem and 

it will be crucial in the definition of avalanche dynamics. A similar equation describing the 

SOC behaviour of extremal models, the gap equation, has been obtained by Bak et al. 0, 

based on phenomenological assumptions. 

The next step is the identification of the scale invariant dynamics for block variables. 

This point can be easily addressed. In fig. ^ we show a coarse graining procedure for the 

extremal dynamics of IP. In the left side of the figure we show some paths leading to cell B (or 

A), each one composed by a set of quenched variables {{u}b)- Each path is characterized by 

the largest variable in the set (e^). The best path leading to cell B is that with the smallest 

efe (saddle point), and will compete with the best path leading to cell A. So, we identify the 

block variables ca and e^ with the saddle point of the best path leading to the corresponding 

cell, say: 

cb = Fb [{u}b\ = min max{ei}b ■ (6) 

b \_ i J 

and analogously for e^- Of course the initial density of the block variables e^ and e^ will 
be rather complex. However, extremal dynamics is not perturbed by a different choice of 
the initial density of the variables. In fact, eq. ^ is invariant under the transformation 
X -^ Jo Po{y)dy which maps the density Po{x) onto the flat one 0. So, we conclude that the 
coarse-grained dynamics is intrinsically scale invariant. 

However, in order to compute the critical IP exponents by real space methods like FST 
we need local growth rules. This corresponds to limit the process inside the FST growth 
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column 1^ , which imphes a truncation of eqs. 0, |^ in order to follow the RTS dynamics of 
only a finite number of variables, instead of those of the whole perimeter. This step must 
be realized carefully, because screening in IP is power law like (we cannot simply ignore 
what happens outside the growth column), and can be realized in the following way. Let 
us consider the asymptotic IP dynamics. Both simulations [^0 and analytical results 0, 
together with eq. ^ show that the fixed point IP dynamics develops itself in scale invariant 
macro events, or avalanches. An avalanche is a spatially and causally connected sequence of 
single growth events. A critical scale invariant avalanche starts with a variable at Pc [H^l- 
By definition, all other variables in the avalanche have values smaller than p^. In order to 
restrict the stochastic process to a limited region of the system we can consider only the 
dynamics inside a single avalanche. In view of the above arguments, the RTS equations for 
this local dynamics are obtained from Eqs. (0, 0) by taking into account only the variables 
which become active after the initiator's growth and by integrating in Eq.(0) only in [0,Pc] 
and not in [0, 1]. This will lead to a GPD for the avalanche dynamics which is defecting, 
i.e. with a normalization less than one, in that it must account of the probability that the 
avalanche stops 0. 

At this point, we use the scale invariant IP local growth rule in the FST scheme, to 
compute the fractal dimension of the infinite percolating cluster [H. This is accomplished 
by computing the transverse nearest-neighbour correlations at a generic scale via suitable 
path integrals. In Table | we show the results of our FST calculations for invasion percolation 
and for some related models. The convergence of the path integrals with respect to the path 
length n is power law like 0], and we need to extrapolate the FST results to n = oo. The 
results of our calculations are in very good agreement both with numerical simulations and 
with known analytical values (Table |). 

We can compute also the exponent of the critical avalanche distribution of IP. The 
distribution of off-critical avalanches in Invasion Percolation {p ^ pc) has been shown to 
have the form 
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D{s;p) = s-^f{\p~p,\sn (7) 

where s is the avalanche size and p is the initiator of an avalanche. In the limit t — > oo 
the system self-organizes into the critical state p = Pc, and the (normalized) avalanche size 
distribution becomes: 

Dis;p,)= ' ^ (8) 

Our calculation scheme develops in the following steps: 1) we evaluate the left hand side of 
Eq. (^; 2) we solve equation (||) for r. For s = 1 we can write: 

D{s = l-pc) = Y. W^^^ ■ P^^\s = l;p,) (9) 

j 

where W^^'' are the weights of the different boundary conditions near the initiator of the 
avalanche and P^^\s = l;Pc) are the stopping probability of the avalanche after one step, 
computed with the local RTS growth rules. 

By Inserting Eq.(|^) in Eq.(^ we get: r ~ 1.5832. Also this analytical result is in 
excellent agreement with numerical simulations [Q, which give r = 1.60. In order to check 
independently our theory, we have performed a new type of computer simulations for the 
avalanche distribution. From RTS scheme and the discussion of the scale invariant local 
dynamics, one notes that /c(to + s), the age of the bond grown during the avalanche at time 
to + s, must satisfy the condition k{tQ + s) < s. So, we can analyze the integer valued signal 
k{t), instead of e(t) P,^] (the value of the smallest variable at time t), in order to estimate 
numerically the exponent r. This alternative method allows us to avoid the problems of 
the numerical approximations that one faces when one analyzes the real signal e{t). We get 
r ~ 1.60 ± 0.02 (fig.^, in very good agreement with our analytical calculation. From the 
knowledge of the exponents Df and r one can recover, via scaling relations [0, all the other 
critical exponents of IP. 

This new theoretical framework, discussed here in relation to the IP problem, can easily 
be extended to other problems of the same extremal nature like, for example, the Bak- 



Sneppen model |13 
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FIGURES 
FIG. 1. Time evolution of the solution of the equation for $((x) (for 2 — d bond IP). ^t{x) 

tends asymptotically to a theta function with discontinuity at pc = 1/2. 

FIG. 2. Renormalization scheme for the extremal dynamics: (a) Dynamics at the smaller 
scale; (b) rescaled dynamics. 

FIG. 3. The distribution of critical avalanches in invasion percolation. The solid line has slope 
-1.60 ±0.02. 
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TABLES 
TABLE L FST results for IP without trapping (Dj), with site trapping (Di), with bond 
trapping {DV), and directed IP {DP^^). The FST results are compared with known analytical 
and simulation values. 



ORDER n 


Df{n) 


D}(n) 


Dy{n) 


DfP{n) 


3 


L7039 


1.6965 


1.7029 


1.6254 


4 


L7941 


1.7378 


1.7825 


1.6626 


5 


1.8228 


1.7506 


1.8066 


1.6924 


6 


L8473 


1.7599 


1.8245 


1.7081 


7 


L8565 


1.7642 


1.8317 


1.7189 


8 


L8645 


1.7678 


1.8372 


1.7250 


9 


L8677 


1.7697 






CXO 


L8879 


1.7812 


1.8544 


1.7444 


analyt. 


§ - 1.895^ 


- 


- 


~ 1.748"^ 


simul. 


~ L89 |] 


~ 1.82 g 


~ 1.86 [| 


- 



^conformal mapping applied to 2d Percolation [12| 
series expansions [11]. 
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